/* Table 2 */
/* be sure mixlogit package is installed */

#delimit;
set mem 30m;
set more off; 

log using Table2.log, t replace;

use formation_new;


mixlogit realg anmax2 pro anti minwin mginvest prevpm gdiv1, 
rand(minor numpar dompar sq median mgodiv1) ln(2) nrep(200) group(case); 

matrix coeffs = e(b);
matrix covmat = e(V); 

keep if case==112;


generate govtstr = subinstr(govt,"govt","",.);

local i = 1;

while `i' < 7 {;

  generate p`i' = substr(govtstr,((`i'*2)-1),2); 
  destring p`i', replace;
  generate party`i' = 0;

  local i = `i' + 1;

};

local k = 1;

while `k' < 7 {;

  replace party`k' = 1 if p1==`k' | p2==`k' | p3==`k' | p4==`k' | 
    p5==`k' | p6==`k';

  local k = `k' + 1;

};

/* rename party dummies */

rename party1 ppr;
rename party2 pvda;
rename party3 d66;
rename party4 vvd;
rename party5 cda;

generate cdavvd = 0;
replace cdavvd = 1 if cda==1 & pvda==0 & vvd==1 & d66==0 & ppr==0;

generate cdapvda = 0;
replace cdapvda = 1 if cda==1 & pvda==1 & vvd==0 & d66==0 & ppr==0;

generate cdad66 = 0;
replace cdad66 = 1 if cda==1 & pvda==0 & vvd==0 & d66==1 & ppr==0;

generate cdpvd66 = 0;
replace cdpvd66 = 1 if cda==1 & pvda==1 & vvd==0 & d66==1 & ppr==0;

generate cdaseat_t = seats if cda==1 & pvda==0 & vvd==0 & d66==0 & ppr==0;
egen cdaseat = max(cdaseat_t);

generate pvdaseat_t = seats if cda==0 & pvda==1 & vvd==0 & d66==0 & ppr==0;
egen pvdaseat = max(pvdaseat);

generate vvdseat_t = seats if cda==0 & pvda==0 & vvd==1 & d66==0 & ppr==0;
egen vvdseat = max(vvdseat);

generate d66seat_t = seats if cda==0 & pvda==0 & vvd==0 & d66==1 & ppr==0;
egen d66seat = max(d66seat);

generate pprseat_t = seats if cda==0 & pvda==0 & vvd==0 & d66==0 & ppr==1;
egen pprseat = max(pprseat);

generate seats2 = seats;

set obs 1000;

generate p1a=.;
generate p2a=.;
generate p3a=.;
generate p4a=.;

generate p1b=.;
generate p2b=.;
generate p3b=.;
generate p4b=.;

generate p1c=.;
generate p2c=.;
generate p3c=.;
generate p4c=.;

generate p1d=.;
generate p2d=.;
generate p3d=.;
generate p4d=.;


drawnorm b1-b19, means(coeffs) cov(covmat) double; 

/* run separate loop for each case to get same draws for eta */

nois _dots 0, title(parametric bootstrap) reps(1000);
forvalues j = 1/1000 {;

  mkmat b1-b19 if _n==`j', matrix(bdraw);
  mixlpred2 pprobs2; 
 
  qui generate p1at = pprobs2 * cdavvd;
  qui generate p2at = pprobs2 * cdapvda;
  qui generate p3at = pprobs2 * cdad66;
  qui generate p4at = pprobs2 * cdpvd66; 


  egen p1att = sum(p1at);
  egen p2att = sum(p2at);
  egen p3att = sum(p3at);
  egen p4att = sum(p4at); 


  qui replace p1a = p1att if _n==`j';
  qui replace p2a = p2att if _n==`j';
  qui replace p3a = p3att if _n==`j';
  qui replace p4a = p4att if _n==`j';


drop   p1at p2at p3at p4at 
       p1att p2att p3att p4att 
       pprobs2;
  
nois _dots `j' 0;
  
};

replace dompar = 0;
replace dompar = 1 if pvda == 1;

replace seats2 = seats;
replace seats2 = seats - cdaseat + pvdaseat if cda == 1;
replace seats2 = seats + cdaseat - pvdaseat if pvda == 1;

replace minor = 0;
replace minor = 1 if seats2 < 0.5;

replace mginvest = minor*invest;

replace minwin = 1;
replace minwin = 0 if seats2 < 0.5;
replace minwin = 0 if (seats2 - pvdaseat)>0.5 & cda==1;
replace minwin = 0 if (seats2 - cdaseat)>0.5 & pvda==1;
replace minwin = 0 if (seats2 - vvdseat)>0.5 & vvd==1;
replace minwin = 0 if (seats2 - d66seat)>0.5 & d66==1;
replace minwin = 0 if (seats2 - pprseat)>0.5 & ppr==1;


nois _dots 0, title(parametric bootstrap) reps(1000);
forvalues j = 1/1000 {;

  mkmat b1-b19 if _n==`j', matrix(bdraw);
  mixlpred2 pprobs2;  

  qui generate p1bt = pprobs2 * cdavvd;
  qui generate p2bt = pprobs2 * cdapvda;
  qui generate p3bt = pprobs2 * cdad66;
  qui generate p4bt = pprobs2 * cdpvd66; 


  egen p1btt = sum(p1bt);
  egen p2btt = sum(p2bt);
  egen p3btt = sum(p3bt);
  egen p4btt = sum(p4bt); 


  qui replace p1b = p1btt if _n==`j';
  qui replace p2b = p2btt if _n==`j';
  qui replace p3b = p3btt if _n==`j';
  qui replace p4b = p4btt if _n==`j';


drop   p1bt p2bt p3bt p4bt 
       p1btt p2btt p3btt p4btt 
       pprobs2;
  
nois _dots `j' 0;
  
};

replace dompar = 0;
replace dompar = 1 if vvd == 1;

replace seats2 = seats;
replace seats2 = seats - cdaseat + vvdseat if cda == 1;
replace seats2 = seats + cdaseat - vvdseat if vvd == 1;

replace minor = 0;
replace minor = 1 if seats2 < 0.5;

replace mginvest = minor*invest;

replace minwin = 1;
replace minwin = 0 if seats2 < 0.5;
replace minwin = 0 if (seats2 - vvdseat)>0.5 & cda==1;
replace minwin = 0 if (seats2 - pvdaseat)>0.5 & pvda==1;
replace minwin = 0 if (seats2 - cdaseat)>0.5 & vvd==1;
replace minwin = 0 if (seats2 - d66seat)>0.5 & d66==1;
replace minwin = 0 if (seats2 - pprseat)>0.5 & ppr==1;


nois _dots 0, title(parametric bootstrap) reps(1000);
forvalues j = 1/1000 {;

  mkmat b1-b19 if _n==`j', matrix(bdraw);
  mixlpred2 pprobs2; 
 
  qui generate p1ct = pprobs2 * cdavvd;
  qui generate p2ct = pprobs2 * cdapvda;
  qui generate p3ct = pprobs2 * cdad66;
  qui generate p4ct = pprobs2 * cdpvd66; 


  egen p1ctt = sum(p1ct);
  egen p2ctt = sum(p2ct);
  egen p3ctt = sum(p3ct);
  egen p4ctt = sum(p4ct); 


  qui replace p1c = p1ctt if _n==`j';
  qui replace p2c = p2ctt if _n==`j';
  qui replace p3c = p3ctt if _n==`j';
  qui replace p4c = p4ctt if _n==`j';


drop   p1ct p2ct p3ct p4ct 
       p1ctt p2ctt p3ctt p4ctt 
       pprobs2;
  
nois _dots `j' 0;
  
};


replace dompar = 0;
replace dompar = 1 if d66 == 1;

replace seats2 = seats;
replace seats2 = seats - cdaseat + d66seat if cda == 1;
replace seats2 = seats + cdaseat - d66seat if d66 == 1;

replace minor = 0;
replace minor = 1 if seats2 < 0.5;

replace mginvest = minor*invest;

replace minwin = 1;
replace minwin = 0 if seats2 < 0.5;
replace minwin = 0 if (seats2 - d66seat)>0.5 & cda==1;
replace minwin = 0 if (seats2 - pvdaseat)>0.5 & pvda==1;
replace minwin = 0 if (seats2 - vvdseat)>0.5 & vvd==1;
replace minwin = 0 if (seats2 - cdaseat)>0.5 & d66==1;
replace minwin = 0 if (seats2 - pprseat)>0.5 & ppr==1;


nois _dots 0, title(parametric bootstrap) reps(1000);
forvalues j = 1/1000 {;

  mkmat b1-b19 if _n==`j', matrix(bdraw);
  mixlpred2 pprobs2;  

  qui generate p1dt = pprobs2 * cdavvd;
  qui generate p2dt = pprobs2 * cdapvda;
  qui generate p3dt = pprobs2 * cdad66;
  qui generate p4dt = pprobs2 * cdpvd66; 


  egen p1dtt = sum(p1dt);
  egen p2dtt = sum(p2dt);
  egen p3dtt = sum(p3dt);
  egen p4dtt = sum(p4dt); 


  qui replace p1d = p1dtt if _n==`j';
  qui replace p2d = p2dtt if _n==`j';
  qui replace p3d = p3dtt if _n==`j';
  qui replace p4d = p4dtt if _n==`j';


drop   p1dt p2dt p3dt p4dt 
       p1dtt p2dtt p3dtt p4dtt 
       pprobs2;
  
nois _dots `j' 0;
  
};

local k = 1;

while `k' < 5 {;

  egen medp`k'a = median(p`k'a);

  egen medp`k'b = median(p`k'b);
  qui generate dp`k'b = p`k'b - p`k'a;
  egen medd`k'b = median(dp`k'b);
  egen lci`k'b = pctile(dp`k'b), p(2.5);
  egen hci`k'b = pctile(dp`k'b), p(97.5);

  egen medp`k'c = median(p`k'c);
  qui generate dp`k'c = p`k'c - p`k'a;
  egen medd`k'c = median(dp`k'c);
  egen lci`k'c = pctile(dp`k'c), p(2.5);
  egen hci`k'c = pctile(dp`k'c), p(97.5);

  egen medp`k'd = median(p`k'd);
  qui generate dp`k'd = p`k'd - p`k'a;
  egen medd`k'd = median(dp`k'd);
  egen lci`k'd = pctile(dp`k'd), p(2.5);
  egen hci`k'd = pctile(dp`k'd), p(97.5);

  egen zmean`k'b = mean(dp`k'b);
  egen zsd`k'b = sd(dp`k'b);
  qui generate z`k'b = zmean`k'b/zsd`k'b;
  qui generate pval`k'b = 2*(1 - normal(abs(z`k'b)));

  egen zmean`k'c = mean(dp`k'c);
  egen zsd`k'c = sd(dp`k'c);
  qui generate z`k'c = zmean`k'c/zsd`k'c;
  qui generate pval`k'c = 2*(1 - normal(abs(z`k'c)));

  egen zmean`k'd = mean(dp`k'd);
  egen zsd`k'd = sd(dp`k'd);
  qui generate z`k'd = zmean`k'd/zsd`k'd;
  qui generate pval`k'd = 2*(1 - normal(abs(z`k'd)));

  local k = `k' + 1;
};

 
summarize p1a p2a p3a p4a 
          p1b p2b p3b p4b 
          p1c p2c p3c p4c
          p1d p2d p3d p4d 
          dp1b dp2b dp3b dp4b 
          dp1c dp2c dp3c dp4c
          dp1d dp2d dp3d dp4d 
          pval1b pval2b pval3b pval4b 
          pval1c pval2c pval3c pval4c 
          pval1d pval2d pval3d pval4d;

summarize medp1a medp2a medp3a medp4a 
          medp1b medp2b medp3b medp4b 
          medp1c medp2c medp3c medp4c
          medp1d medp2d medp3d medp4d 
          medd1b medd2b medd3b medd4b 
          medd1c medd2c medd3c medd4c
          medd1d medd2d medd3d medd4d 
          lci1b hci1b lci2b hci2b 
          lci3b hci3b lci4b hci4b 
          lci1c hci1c lci2c hci2c 
          lci3c hci3c lci4c hci4c 
          lci1d hci1d lci2d hci2d 
          lci3d hci3d lci4d hci4d;

